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Abstract 

Background The cyclic AMP- Protein Kinase A (cAMP-PKA) pathway is an evolutionarily conserved signal transduction 
mechanism that regulates cellular growth and differentiation in animals and fungi. We present a mathematical model 
that recapitulates the short-term and long-term dynamics of this pathway in the budding yeast, Socchoromyces 
cerevisiae. Our model is aimed at recapitulating the dynamics of cAMP signaling for wild-type cells as well as single 
{pdelA and pde2A) and double {pdel A pde2A) phosphodiesterase mutants. 

Results Our model focuses on PKA-mediated negative feedback on the activity of phosphodiesterases and the Ras 
branch of the cAMP-PKA pathway. We show that both of these types of negative feedback are required to reproduce 
the wild-type signaling behavior that occurs on both short and long time scales, as well as the the observed responses 
of phosphodiesterase mutants. A novel feature of our model is that, for a wide range of parameters, it predicts that 
intracellular cAMP concentrations should exhibit decaying oscillatory dynamics in their approach to steady state 
following glucose stimulation. Experimental measurements of cAMP levels in two genetic backgrounds of 5. cerevisiae 
confirmed the presence of decaying cAMP oscillations as predicted by the model. 

Conclusions Our model of the cAMP-PKA pathway provides new insights into how yeast respond to alterations in 
their nutrient environment. Because the model has both predictive and explanatory power it will serve as a 
foundation for future mathematical and experimental studies of this important signaling network. 
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Background 

In eukaryotic cells, the cyclic adenosine monophosphate 
(cAMP) - Protein Kinase A (PKA) pathway plays a cen- 
tral role in mediating diverse biological responses such as 
growth, development, and cell differentiation [1,2]. Stim- 
uli such as hormones, neurotransmitters, nutrients, and 
physiological stress agents trigger signaling cascades that 
collectively mount a response through cAMP-mediated 
PKA signaling [3-9]. In the budding yeast, Saccharomyces 
cerevisiae, activation of PKA in response to essential 
nutrients and fermentable carbon sources is directed by 
intracellular levels of cAMP [10]. Synthesis of cAMP 
from ATP is catalyzed by the enzyme adenylate cyclase 
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(Cyrl) and governed by two different G-protein systems 
(Figure 1) [10-12]. Ras2, a small GTP-binding protein, in 
its G TP-bound (active) state stimulates adenylate cyclase, 
causing a rapid increase in intracellular cAMP levels 
(Figure 1, top right) [11]. Ras2 activity is positively regu- 
lated by the guanine nucleotide exchange factors (GEFs) 
Cdc25 and Sdc25 as a function of intracellular glucose lev- 
els [13-20]. Conversely, GTPase activating proteins, Iral 
and Ira2, down regulate Ras2 activity by stimulating the 
hydrolysis of GTP [21,22]. In parallel to Ras, a second 
G-protein pathway, involving the proteins Gprl, Gpa2, 
and Rgs2, responds to extracellular levels of glucose and 
increases adenylate cyclase activity (Figure 1, top left). 
Gprl is a membrane bound G-protein coupled receptor 
that activates the Ga protein Gpa2 in response to extracel- 
lular glucose levels. The GTP bound form of Gpa2 in turn 
stimulates adenylate cyclase activity [12,20,23,24], cAMP 
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Figure 1 Key interactions of the cAMP-PKA pathway in yeast. See the main text and Additional file 1 : Supplementary Model for descriptions of 
each reaction. Arabic numerals refer to equations in the Additional file 1 : Supplementary Model that describe each reaction. The key interactions in 
the model presented here are depicted with red lines. 



synthesis through this pathway is negatively controlled by 
Rgs2, the GTPase activating protein of Gpa2 [25,26]. An 
increase in cellular cAMP levels upon glucose induction 
activates cAMP-dependent PKA. PKA is a holoenzyme 
that includes both regulatory (Bcyl) and catalytic sub- 
units (Tpkl, Tpk2 and Tpk3) [27,28]. cAMP binding to 
the regulatory subunits leads to the release of the catalytic 
PKA subunits ([29]; Figure 1, bottom) which are then free 
to interact with downstream targets such as metabolic 
enzymes, transcription factors and other kinases [30-33]. 

Downstream of the adenylate cyclase, cAMP levels are 
modified by the action of phosphodiesterases, enzymes 
that catalyze the conversion of cAMP to AMP, thus pre- 
venting cAMP accumulation in the cell [34,35]. In yeast 
there are two phosphodiesterases, the low-affinity phos- 
phodiesterase Pdel and the higher affinity Pde2. Pde2 
is a class II phosphodiesterase and shares homology 



with other eukaryotic phosphodiesterases [36], while Pdel 
homologues have been identified in only a small number 
of eukaryotes [35]. The high-affinity Pde2 is thought to 
play the key role in maintaining cAMP at steady state lev- 
els [37]. The role of Pdel is less well understood, though 
genetic studies in both S. cerevisiae and the distantly 
related yeast Candida albicans suggests that Pdel is a tar- 
get of PKA, and that PKA-mediated feedback on Pdel is 
an important component of maintaining tight regulation 
ofcAMP levels [35,38]. 

In this study we present a mathematical model of the 
cAMP-PKA pathway in Saccharomyces cerevisiae. The 
primary goal of this model is to explain the observed 
response of wild- type cells and PDE mutants (pdel A, 
pde2A, and pdelApde2A) to glucose stimulation. Wild- 
type cells exhibit a rapid, transient increase in cAMP 
levels followed by a quick return to a new steady state; 
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this response typically occurs on a time scale of 60-90 
seconds [35]. pde2A mutants have approximately wild- 
type dynamics, while pdelA mutants exhibit elevated 
peak cAMP levels and a slow return to a steady state 
that is higher than that of wild- type cells. Surprisingly, 
the double phosphodiesterase mutant, pdelApde2A, 
shows no response to glucose stimulus, and main- 
tains cAMP concentrations at essentially a pre-stimulus 
steady state [35]. The question we posed was whether 
it was possible to mathematically reconcile all of these 
observed dynamics with currently known genetic and 
biochemical interactions. Our model focuses on PKA- 
mediated negative feedback interactions, and suggests 
that feedback on both phosphodiesterases and on the 
Ras branch of the pathway are important for reproduc- 
ing the dynamics observed in wild-type cells and PDE 
mutants. 

While a number of studies have proposed mathematical 
models for the cAMP-PKA network in yeast [39-41], none 
of these previous models is capable of recapitulating the 
observed dynamics of both wild-type cells and the double 
phosphodiesterase mutant (pdelA pde2A). Our proposed 
model is thus unique in this aspect. Another notable fea- 
ture of our model is that, for a wide range of parameters, 
it predicts the presence of decaying cAMP oscillations fol- 
lowing glucose stimulus. Oscillatory behaviors were not 
an explicit goal of our initial modeling effort, and thus we 
treated this as a novel model prediction. We tested this 
prediction by measuring the cAMP response to glucose 
stimulus in two genetically diverse yeast strains (S288c 
and £ 1278b). Our experiments reveal variation in oscil- 
latory behaviors between these genetic backgrounds, and 
we demonstrate that our model, in turn, is capable of 
reproducing the diversity of observed oscillatory patterns. 
Finally, our core model can also be easily extended to gen- 
erate sustained cAMP oscillations that occur on longer 
time scales (tens of minutes; [42]). The results presented 
here thus demonstrate both the descriptive and predictive 
power of the model we have developed, and suggest that 
this model can serve as a basis for future experimental 
and modeling studies of this important eukaryotic signal 
transduction pathway. 

Summary of motivating experiments 

There is a large body of experimental work that provides 
information on genetic and biochemical interactions rele- 
vant to the cAMP-PKA pathway. The strategy we adopted 
in this study was to use this body of work as a basis 
for constraining interactions in our model. In particu- 
lar our model was developed to understand two specific 
experimental results: the short-term dynamics of cAMP 
following glucose stimulation reported by Ma et al. [35]; 
and the long-term dynamics of cAMP under various stress 
levels reported by Garmendia-Torres et al. [42]. 



Short-term cAMP dynamics 

Ma et al. [35] investigated the short-term behavior of 
cAMP signaling following glucose stimulation in cells 
that had previously been grown under carbon-source 
limitation. These data therefore represents the dynam- 
ics of cAMP signaling as cells transition from starved 
to fed states, on a time scale of seconds to minutes. 
Ma et al.s analysis focused particularly on PKA-mediated 
negative feedback involving the phosphodiesterases Pdel 
and Pde2. Figure 2 summarizes the five dynamic pat- 
terns we consider here, based on Ma et al.s experiments 
involving wild-type cells and various phosphodiesterase 
mutants (compare to Ma et al. Figures 2a and 3a). These 
cases are: 

Case 1: Wild-type (wt; blue line) - the concentration of 
cAMP in wt cells rapidly increases following 
glucose stimulation. After reaching a peak level, 
the concentration of cAMP declines to a new 
steady state that is higher than its initial 
concentration. We refer to this as the wild-type 
transient response. 

Case 2: Pdel knockout (pdel A; red line) - the 

concentration of cAMP exhibits a much larger 
and longer transient response, about four times 
as large as the transient response exhibited by 
wild-type. 

Case 3: Pde2 knockout (pde2 A; green line) - the 
concentration of cAMP exhibits similar 
dynamics to wild-type. 

Case 4: Pdel phosphorylation mutant (pdel ala152 ; pink 
line) - the cAMP transient response is about 
twice as large as that exhibited by wild-type. 

Case 5: Double PDE knockout (pdel A pde2A; black 
line) - glucose stimulation has no effect on the 
concentration of cAMP present in the cytosol. 
The initial value in this case is much larger then 
the wild-type initial value. 

These results are surprising in two ways: 1) the dynam- 
ics observed in the pde2A mutant are almost the same as 
those observed in wild-type cells, while the response in 
the pdel A mutant is drastically different, despite the fact 
that Pde2 is thought to be the higher affinity phospho- 
diesterase; and 2) in the double mutant (pdelApde2A) 
glucose stimulation has no effect on the concentration 
of cAMP, indicating a synergistic interaction between the 
PDE null-alleles. 

Long-term cAMP dynamics 

Garmendia-Torres et al. [42] focused on the long-term 
behavior of cAMP signaling in yeast exposed to stress 
conditions. They predicted sustained oscillations, with a 
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Figure 2 cAMP dynamics following glucose stimulus as described by Ma et al. wild type (blue); pdel A (red), pde2A (green), 
p£ /e7 ala152 (pink) ,pde1Apde2A (black) . 



period of 5-10 minutes, in the concentration of cAMP 
at intermediate stress levels. This prediction was moti- 
vated by oscillations observed in the nucleocytoplasmic 
shuttling of the transcriptional activator Msn2. Msn2 is 
negatively regulated by the cAMP-PKA pathway through 
phosphorylation of Msn2 by PKA; this phosphoryla- 
tion promotes export of Msn2 from the nucleus [43]. 
Garmendia-Torres et al. predicted that the nucleocyto- 
plasmic shuttling of Msn2 that they observed was driven 
by oscillations in the concentration of cAMP, and they 
constructed a mathematical model to explore the pro- 
posed dynamics. In their model, high levels of stress 
drive the concentration of cAMP to a low steady state, 
and at low stress levels the cAMP levels go to a high 
steady state. Only under intermediate levels of stress did 
cAMP oscillations occur. This oscillatory behavior in their 
model was driven by PKA-mediated negative feedback on 
RasGTP. 

Mechanisms of PKA-mediated negative feedback on cAMP 
signaling 

Both sets of observations described above suggest that 
PKA plays an important role in attenuating cAMP sig- 
naling via negative feedback interactions. PKA-mediated 
negative feedback on cAMP signaling has been proposed 
to occur at two levels - via enhancement of phosphodi- 
esterase activity and via attenuation of Ras signaling. 

PKA interactions with Pdel are supported by experi- 
ments by Ma et al. that showed that mutagenesis of a 
putative PKA phosphorylation site in Pdel (pdel ala152 ) 



causes a dramatic increase in cAMP accumulation 
(Case 4, Figure 2). There is less direct evidence for PKA 
interactions with Pde2 but other class II phosphodi- 
esterases are known to be targets of PKA [36]. Hu et al. 
[44] showed that hyperactive PKA activity leads to ele- 
vated Pde2 levels, and that this resulted from increased 
protein stability of Pde2. 

The observation that cAMP signaling did not signifi- 
cantly change in the pde2A mutant might seem to suggest 
that Pde2 plays no role in regulating cAMP immediately 
following glucose stimulus. However, the cAMP transient 
in the pdel A mutant, and the synergistic effect observed 
in the pdel A pde2A double mutant, argue for a role for 
Pde2, though with a slower activation rate as compared 
to Pdel. However, while we favor the hypothesis of PKA 
feedback on both PDEs, we can not rule out the possibility 
that the transient cAMP dynamics observed in the pdel A 
mutant might involve other feedback mechanisms. 

With respect to Ras signaling, one or more elements 
of this pathway may be targets of PKA-mediated feed- 
back. For example, Iral and Ira2 both have predicted 
PKA phosphorylation sites [12,45]. Similarly, PKA hyper- 
phosphorylates Cdc25, causing its dissociation from 
the adenylate cyclase/Ras2 complex. Cdc25 may also 
directly bind to the adenylate cyclase to facilitate mem- 
brane anchoring [46-48]. Re-localization of Cdc25 to the 
cytoplasm attenuates the glucose response of adenylate 
cyclase, inhibiting synthesis of cAMP [49,50]. More- 
over, membrane dissociated Cdc25 can no longer acti- 
vate Ras2, triggering an increase in Ras-GDP levels [4]. 
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Ma et al suggested that the likely explanation for the pat- 
tern they observed in the pdelApde2A double mutant 
(no change in cAMP levels after glucose stimulation) 
resulted from constitutively high PKA feedback inhibition 
in the absence of the phosphodiesterases. 

In summary, both studies described above suggest that 
PKA down regulates cAMP signaling by increasing the 
rate of PDE mediated decay and/or by decreasing the 
rate of cAMP production via interactions with the Ras 
effectors. Thus a key goal of our modeling effort was 
to explore PKA feedback on phosphodiesterase activity 
and the Ras pathway. In our model PKA interacts with 
both PDEs to increase their activity, though we assume 
PKA has a higher affinity for Pdel than for Pde2. Sim- 
ilarly, we model negative feedback on Ras-GTP through 
PKA phosphorylation of the Ras-GTPases Iral and Ira2. 
These interactions are highlighted by the red arrows in 
Figure 1. 

Model 

To fully describe the key reactions of the cAMP-PKA 
pathway as depicted in Figure 1 requires a system of 
nine time-dependent differential equations to model the 
concentrations of: 1) Gpa2GTP (active Gpa2); 2) Cdc25; 
3) Ras-GTP-ases Iral and Ira2; 4) Ras-GTP (active Ras2); 
5) adenylate cyclase activity; 6) activated Pdel; 7) acti- 
vated Pde2; 8) cAMP; and 9) active PKA (free catalytic 
subunits). The full nine-dimensional system is described 
in detail in Additional file 1: Supplementary model 
(Equations (S1)-(S9)). 

For the purposes of exposition, we first describe the 
motivation behind the equations that represent PKA 
mediated feeedback. We then describe a simplification of 
the model to a four-dimensional system that replicates the 
key dynamical behaviors of cAMP signaling as described 
above. 

Key reactions for PKA-mediated feedback 

PKA feedback (red arrows in Figure 1) in our model 
occurs in three ways: (i) PKA phosphorylates Pdel, 
enhancing Pdel breakdown of cAMP; (ii) PKA phospho- 
rylates Pde2, enhancing Pde2 breakdown of cAMP; and 
(iii) PKA phosphorylates the Ras-GTP-ases Iral and Ira2. 
In Cases (i) and (ii), PKA inhibits the concentration of 
cAMP by enhancing the breakdown of cAMP. In Case 
(iii) PKA inhibits the production of cAMP by inactivating 
Ras-GTP. 

We assume that there are much greater concentra- 
tions of Pdel, Pde2, Iral and Ira2 in the cytoplasm than 
of PKA. We also assume that Pdel, Pde2 and Iral/2 
compete with each other for activation by PKA. We fur- 
ther assume that Iral and Ira2 can be simultaneously acti- 
vated by PKA, thus we treat them as a single variable, Ira. 



These assumptions are taken in order to understand the 
Ma et al. cases. The three competitive reactions are spec- 
ified as follows: 

[ Pdel] + [ PKA] ^[Pdel :: PKA] %[ Pdel*] +[ PKA] 

k-i 

[ Pde2] + [ PKA] ^ [ Pde2 :: PKA] -5 [ Pde2*] +[ PKA] 

k-2 

[ Ira] + [ PKA] h [ Ira :: PKA] h [ Ira*] +[ PKA] 

k-z 

Here the asterisk (*) indicates the activated or phos- 
phorylated form of the enzyme. Here and below, PKA 
refers the catalytic subunits (Tpkl, Tpk2, or Tpk3). 
Equations (l)-(3), below, model the velocities of these 
reactions using an extension of Michaelis-Menten kinetics 
[51] that allows for competition between Pdel, Pde2 and 
Ira for PKA. 



Pdel* 



7?^ [Pdel] [PKA] 



dt 



4-[Pde2*]: 
dt 



d - 

dt 
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i>+[Pdel] +^[Pde2] +T^[Ira] 

#p 2 [Pde2] [ PKA] 
r p2 +[ Pde2] +{g [ Pdel] + ^[ Ira] 

. y? z [Ira][PKA] 
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(1) 



(2) 



(3) 



In order for the model to replicate the dynamics 
observed by Ma et al. [35] (Figure 2) we impose the follo- 
wing conditions: 

Condition (a) The following inequalities must hold: 
I>[Pde2] 



I>[Pdel] 
I> 2 [Ira] 



«1 



T,[Pde2] <<C1 

Condition (b) In comparing analogous reactions of 

Pdel and Pde2, the reactions of Pde2 are 
uniformly slower. 

Condition (c) PKA rapidly phosphorylates Ira. 

By invoking Condition (a), the effect of Pde2 is negli- 
gible when Pdel is active. This assumption is justified by 
the observation that the dynamics of cAMP signaling fol- 
lowing glucose stimulus are essentially unchanged in the 
pde2A mutant. By invoking Condition (b) we can model 
both Pdel and Pde2 with a single equation (for details 
and derivation see the Additional file 1: Supplementary 
Model). Condition (b) is justified by the large transient 
cAMP signal observed in the pdel A mutant. Neither 
Condition (a) nor (b) are required to reproduce the 
behavior of the wild-type or mutant strains, but are 
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mathematically convenient for exploring the model ana- 
lytically as described below. 

Condition (c) reflects negative feedback on the Ras 
branch of the pathway. We impose Condition (c) so that 
even when PKA's phosphorylation of Ira is slowed by the 
presence of either Pdel or Pde2, by Condition (a), the 
effect of PKA on Ira can not be ignored. In the double 
mutant case (pdelApde2A) the loss of competitive inhi- 
bition of PKA greatly enhances PKAs effect on Ira. Thus 
in the double mutant case PKA has a strong negative effect 
on Ras-GTP. Several studies have provided experimental 
support for a role of the Ira proteins in negative feedback 
on cAMP signaling [21,22] as well as evidence that this 
feedback is regulated by PKA [52] . 

Steady state assumptions 

The model may be further simplified by assuming that 
the following four reactions are fast and hence proceed to 
steady state. 

1. Gpa2-GDP ^ Gpa2-GTP (reaction 1 in Figure 1). 

2. Activation/inactivation of Cdc25 (reaction 2 in 
Figure 1). 

3. Activation/inactivation of adenylate cyclase (reaction 
5 in Figure 1). 

4. cAMP+[Bcyl::PKA] ^ [cAMP::Bcyl] +PKA 
(reaction 9 in Figure 1). 

Again, these steady state assumptions are made for 
mathematical simplicity, and are not a requirement to 
mathematically reproduce wild-type or mutant signaling 
dynamics. We will only focus on the last of these steady 
state assumptions, for the rest we refer the reader to the 
Additional file 1: Supplementary Model. We assume that 
PKA activates rapidly, that is we assume that the concen- 
tration of active PKA can be treated as being at steady 
state. This assumption was taken for mathematical sim- 
plicity of the model, it is not a necessary condition to 
numerically replicate the dynamics observed by either Ma 
et al. or Garmendia-Torres et al. In our model PKA feed- 
back plays a central role, thus we will examine this steady 
state condition in detail. In the cAMP-PKA pathway PKA 
is activated by cAMP in the following manner: four cAMP 
molecules bind to two regulatory subunits (Bcyl) and 
release two catalytic subunits, creating the active form 
of PKA. 

4[cAMP] + [2Bcyl ::PKA] ^[4cAMP:: 2Bcyl] + [ 2PKA*] 

k f 

When modeling this reaction we can simplify this in 
two ways: 1) by assuming that the concentration of active 
PKA does not approach its maximum; and 2) by approx- 
imating the 4 to 2 ratio by a 2 to 1 ratio. Applying this 
simplification we can model the change in PKA activity 



as proportional to cAMP squared (x 2 ); thus in Equations 
(4a) -(4d) (below) x 2 represents PKA feedback (for full 
details see the Additional file 1: Supplementary Model). 

Simplified model 

Combining feedback Conditions (a), (b), and (c) with 
the four steady state assumptions, and by modeling both 
forms of PDE with a single case dependent variable, p, we 
reduce the nine-dimensional model shown in Figure 1 to 
the simpler four-dimensional model depicted in Figure 3. 
The model in Figure 3 is described by the following system 
of equations: 



Ras : 
Ira : 
Pde: 
cAMP : 



dr 
dt 
dz 

dt 
dp 

dt 
dx 

dt 



A(l - r) Bzr 
ri + 1 - r ri + r 

= N(x 2 - z) 



= M(x 2 -p) 



C+Gr- D 0 x - 



Dpx 

r +x 



(4a) 
(4b) 
(4c) 
(4d) 



In Equations (4a) -(4d) and in Figure 3 the concentration 
of Ras-GTP is represented by variable r; the concentration 
of Ras -GDP is represented by (1 — r); the concentrations 
of Iral and Ira2 are represented by z; the variable p rep- 
resents the concentration of Pdel in Cases 1, 3, and 4 (wt, 
pde2A, and pdel ala152 ) and the concentration of Pde2 in 
Case 2 (pdel A); the concentration of cAMP is represented 
by x; as stated above x 2 represents the concentration of 
active PKA; and time is given by variable t All vari- 
ables have been nondimensionalized. The constants in 
Equations (4a)-(4d) are nondimensional composites of 
dimensional parameters of the full system, derived in the 
Additional file 1: Supplementary Model. Intuitively, A is 
the activation rate of Ras-GTP catalyzed by Cdc25; B is the 
inactivation rate of Ras-GTP catalyzed by Iral and Ira2; 
Ti is the affinity of both Cdc25 and Iral/2 for Ras-GTP; 
N is the reaction rate of Iral and Ira2; M is the quantita- 
tive expression of Condition (b); C is the production rate 
of cAMP due to the basal activity of adenylate cyclase; G 
accounts for the glucose, normalized so that G = 1 after a 
glucose stimulus is applied; Do represents a "basal" decay 
rate of cAMP in the absence of activated Pde; D represents 
enhanced decay due to PKA feedback on Pde modified by 
T, the affinity of activated Pde for cAMP. 

Results 

In this section we examine analytical, numerical, and 
experimental results motivated by our model of the 
cAMP-PKA pathway. We show that the simplified version 
of the model, given by Equations (4a) -(4d), can adequately 
replicate the short-term dynamics of cAMP, reported by 
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Ma et al [35] for wild-type cells and all three phospho- 
diesterase mutants (pdelA, pde2A, pdelApde2A). For 
the wild-type case we analyze the behavior of the sys- 
tem as it approaches steady state and find that for a 
wide range of parameters the model predicts that cAMP 
levels should exhibit decaying oscillations. Since oscilla- 
tory behavior was not an explict input into our mod- 
eling effort, we considered this a novel prediction, and 
undertook an experimental validation of the model using 
time series measurements of cAMP signaling for two 
diploid yeast strains. As predicted by the mathematical 
model we observed decaying cAMP oscillations follow- 
ing glucose stimulation. We then extended the model to 
consider longer-time scale, sustained cAMP oscillations 
(5-10 minute periods) in response to nutrient stress, as 
predicted by Garmendia-Torres et al. [42]. We show that 
the core model can also reproduce sustained oscillations 
on this longer time scale, and we discuss the model param- 
eters and the corresponding biochemical interactions, 
that are required to generate such sustained oscillations. 

Short-term dynamics 

Our first task is to demonstrate that Equations (4a)- 
(4d) can adequately replicate the dynamics of cAMP 
signaling immediately following glucose stimulation, for 
both wild-type cells and mutants lacking either one or 
both phosphodiesterases. Figure 4 shows the modeled 
concentration of cAMP as a function of time, follow- 
ing glucose stimulation of wild-type and PDE mutants. 



This outcome was generated using the parameter choices 
given in Tables 1 and 2. The dynamical patterns gen- 
erated by the model are a good match to those illus- 
trated in Figure 2. The dimensional concentrations in 
the figure were obtained by multiplying x in Equations 
(4a)-(4d), the dimensionless variable representing cAMP 
concentrations, by 24.95 fmol • 10 -6 cells; and multiplying 
time, given by dimensionless variable t, by .038 minutes. 
These scale factors emerge from the discussion in the 
Additional file 1: Supplementary Model. Initial conditions 
were the steady-state concentrations that occur under 
glucose starved conditions. 

Note that the numbers in columns of Table 1 that repre- 
sent the wild-type and pde2A backgrounds are identical. 
Because of Condition (a) above, when Pdel is present, 
the effect of Pde2 is negligible. Thus, the values in these 
columns of the table characterize properties of Pdel. Sim- 
ilarly, the values in the pdel A column characterize prop- 
erties of Pde2. In the case of pdel ala152 both the activation 
and inactivation rate of Pdel decrease, but the affinity of 
Pdel for cAMP (r) does not change. The value of B is 
increased in the double PDE knockout because, by Con- 
dition (a), when both forms of PDE are eliminated the 
effect of PKA on Ira 1/2 is greatly enhanced. As men- 
tioned above, the dimensionless parameters in Equations 
(4a) -(4d) are composites of dimensional parameters. In 
choosing dimensional parameters we used values given by 
Garmendia-Torres et al. [42] when available. Moreover, 
the parameters determined by fitting the model also lie in 
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Figure 4 Numerical simulations under the simplified model Equations (4a)-(4d), fit to the data of Ma et al (Parameters as given in 
Tables 1 and 2). The corresponding genotypes are as in Figure 2. The wt and pde2A cases overlap each other in this figure. 



biologically appropriate ranges. We refer to the Additional 
file 1: Supplementary Model for the details of the fitting 
process. 

An intuitive explanation for the behavior of wild-type 
cells is that the transient cAMP peak is due to the 
delay between the increase of adenylate cyclase activity 
and PKA mediated enhancement of phosphodiesterase 
activity. Activation of adenylate cyclase activity increases 
cAMP concentrations, leading to the release of PKA cat- 
alytic subunits. PKA in turn activates the phosphodi- 
esterases, which enhances the decay of cAMP and brings 
the concentration of cAMP down to its steady state level. 
In wild type cells, we propose that Pdel is the primary 
effector of PKA feedback. In the pdel A and pdel ala152 
cases the slower activation of Pde2 accounts for the 
greater transient observed in these cases. 



Oscillatory approaches to steady state 

For a wide range of parameters the model described by 
Equations (4a) -(4d) suggests that cAMP concentration 
should exhibit decaying oscillations as it approaches its 
steady state value following glucose stimulus. Here we 
describe the conditions required to generate this behavior. 

Decaying cAMP oscillations via feedback on PDEs 

PKA feedback through the PDEs alone is sufficient to 
generate decaying cAMP oscillations following glucose 
stimulus. To see this we examine the case when Ras is 
either completely in the Ras-GTP or Ras-GDP state. When 
A is sufficiently large, or sufficiently small, we can reduce 
the model by making the approximation that the concen- 
tration of Ras stays completely either in the G TP-bound 
(r ss & 1), or GDP-bound form (r ss « 0). This reduces 



Table 1 Parameter values that vary when fitting 
Equations (4a)-(4d) to the Ma et al. data 



Parameter 



wt 



pc/elA pde2A pde\ 



alal 52 



pde'\A 
pdel A 



Here M is the quantitative expression of Condition (b) [see text]; D represents 
enhanced decay due to PKA feedback on Pde modified by r, the affinity of 
activated Pde for cAMP; B is the inactivation rate of Ras-GTP due to Ira 1/2. Also 
see the text for a discussion of the great overlap between columns 2 and 4. 



Table 2 Parameter values that are consistent when fitting 
Equations (4a)-(4d) to the Ma et al. data 



Parameter 



Value 



M 


0.01 


0.0005 


0.01 


0.0025 




Ti 


D 


1 


0.26 


1 


0.54 




N 


r 


33.6 


16.8 


33.6 


33.6 




C 


B 


.0051 


.0051 


.0051 


.0051 


.51 


Do 



1.45 
0.0004 
0.032 
0.044 
0.013 



Here A is the activation rate of Ras-GTP catalyzed by Cdc25; r"i is the affinity of 
both Cdc25 and Ira 1/2 for Ras-GTP; N is the reaction rate of lral/2; C is the 
production rate of cAMP due to the basal activity of adenylate cyclase; D 0 
represents the decay rate of cAMP in the absence of activated PDEs. 
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our model, Equations (4a)- (4d), to a two-by-two sys- 
tem, representing change in concentration of Pde (p) and 
cAMP (x): 



^=M(x 2 -p) 

dx Dpx 
— = C 0 -D 0 x- ——, 
at r + x 



(5a) 



(5b) 



where Q = 1 + C when Ras-GTP is the dominant Ras 
form (r ss ~ 1) or Q = C when Ras -GDP dominates 
(r ss & 0). In these cases we can find an explicit expression 
for when the concentration of cAMP exhibits decaying 
oscillation as it approaches steady state. 

Claim L Solutions of(Sa) and (Sb) exhibit an oscillatory 
approach to steady-state if and only if 



^ TDx 2 ss 
M-D 0 - 



Dxl 
< 8M- ss 



\ (r h~ %ss) ) r + x ss 

where x ss is the equilibrium solution of (5a, Sb); i.e., 

Dx i 

Cq — DqX ss — 



r + x sl 



■ 0. 



(6) 



(7) 



Condition (6) is derived by computing the Jacobian 
matrix of Equations (5a) and (5b) at the equilibrium point 
and determining when the eigenvalues have a nonzero 
imaginary part. Because of the implicit definition of x ss , 
it is rather difficult to apply directly. Therefore we refer 
to Figures 5a and 5b, in which the overlapping shaded 
regions show the ranges of Do and D for which the 
inequality (6) is satisfied, for several values of T and for 
M assumed equal to .01. Figure 5a shows the parameter 
range for Co = 1 + C = 1.044, the case when r ss ~ 1, and 
Figure 5b shows the parameter range for Co = C = .044, 
the case when r ss & 0. Observe that the value of V has 
a large effect on when oscillation occurs: when V is small 
(i.e. the PDEs have high affinity for cAMP), oscillations 
occur over a large parameter range, but as T increases, the 
parameter range decreases significantly. Finally, because 
of the factor of M on the right-hand side of (6), it is dif- 
ficult to satisfy the inequality if M is very small; thus 
decaying oscillations are less likely for the pdelA and 
pdel ala152 mutants. 

General conditions for decaying oscillations 

We now seek to understand the parameters that cause 
decaying cAMP oscillations when Ras activation is free 
to vary. To do so we numerically examine when the Jaco- 
bian matrix determined by the four-by-four system has a 
complex conjugate pair of eigenvalues with negative real 
parts. We summarize the results here. We find that as 
parameters Ti, N and the ratio ^ increase we are more 
likely to observe decaying oscillations. We also find that 
as the parameters that enhance the decay of cAMP are 



increased, that is as parameters Dq and D are increased 
or as M decreases decaying oscillations are less likely. 
For a more detailed analysis we refer the reader to the 
Additional file 1: Supplementary Model. 

Experimental validation of predicted cAMP oscillations 

Our model makes the novel prediction that for many 
choices of parameter values, cAMP levels should exhibit 
decaying oscillations towards a steady-state following glu- 
cose stimulus. In order to test this hypothesis we mea- 
sured the dynamics of the cAMP response in diploid cells 
of two strains of S. cerevisiae - S288c and E 1278b. S288c 
is considered the standard "reference" genome for yeast 
studies while E 1278b is commonly used for studies of 
developmental pathways in yeast [53,54]. 

For each strain we monitored cAMP levels for eight to 
twelve minutes following a glucose stimulus. Two typi- 
cal experimental time series are illustrated in Figure 6. 
There are cAMP oscillations in both genetic backgrounds 
but the quantitative features of the oscillations appear to 
be strain dependent. E 1278b exhibits a classical form of 
decaying oscillations characterized by a large cAMP peak 
immediately after the stimulus, and dampening oscilla- 
tions towards a new steady state. In the S288c background, 
by contrast, the oscillations appear to be delayed, occur 
with a lower amplitude and do not decay as rapidly. The 
E 1278b strain has a cAMP peak approximately twice as 
large as that of S288c. 

There are various ways in which we could fit the dynam- 
ics reported here. We could assume that Ras-GTP is 
saturated, as in Ma et al. Cases 1-4 and seek to fit only 
the two-by-two system, given by Equations (5a) and (5b), 
to this data. Alternatively we could fit all four equations, 
(4a) -(4d), to this data. Since this experiment was done in 
a similar way to the Ma et al. [35] experiment, we seek to 
fit the data to parameters relating to only Pde feedback, 
that is parameters D, T, and M. The data is taken in incre- 
ments of fifteen seconds for the first two minutes, thirty 
seconds for the next five minutes and every minute for the 
last five minutes. We used a least squares approach to fit 
the parameters of the mathematical model to the observed 
data. Because the density of samples is higher in the initial 
two minutes, the fitting more closely resembles the data 
over this interval. 

We can approximate the dynamics observed for both 
E 1278b and S288c using the model described by 
Equations (4a) -(4d) (Figure 6A). Table 3 shows the param- 
eter values used to fit the oscillations in E 1278b and S288c 
compared to the parameter values used to fit Ma et al. 
wild type. For both cases the model requires Pdes affin- 
ity for cAMP, ^, to be much greater than Pde s affinity for 
cAMP in the Ma et al. case. In our fitting we note that the 
most significant difference between the two strains is in 
the decay of cAMP with respect to Pde (parameter D). 
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Figure 5 Solutions to the inequality given in Equation (6) for two different values of the model parameter C. A) C 0 = 1 044; B) C 0 = 044 The 

overlapping shaded regions of these figures illustrate the parameter ranges for D and Do where the model predicts the possibility of cAMP oscillation. 
The black region is r = 33.6, the r value for Ma et al. wt; the blue region is r = 3.36; and the red region is r = .336. In both cases, M = .01 . 



PDE mutants in the £ 1278b background 

The £ 1278b wild-type strain exhibited a considerably 
larger cAMP peak and more pronounced oscillations than 
S288c (Figure 6), or the the W303 backround used by Ma 
et al We therefore undertook additional experiments to 
examine the behavior of the pdelA and pde2A knockout 
mutants in £ 1278b. Typical experimental time series for 
the PDE mutants in the £ 1278b background are shown in 
Figure 7. These experiments qualitatively agree with the 
dynamics reported by Ma et al for the W303 background: 



1) the pde2A mutant has very similar cAMP dynamics 
to the wild-type background; and 2) the pdelA mutant 
reaches a much higher cAMP maximum. One notable dif- 
ference between the W303 and £ 1278b pdelA mutants 
is that cAMP levels in the E 1278b background reach a 
steady-state at or near their maximum, while the W303 
pdelA mutant shows a more transient response with 
cAMP levels decaying towards a steady state significantly 
lower than the maximum. Because of this observed differ- 
ence in the cAMP signaling behavior, when fitting model 
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Figure 6 Variation in cAMP signaling between strain backgrounds. A) Experimental time series of cAMP concentration following glucose 
stimulus for S288c (red) and £ 1 278b (blue), x's indicate measured cAMP values, the solid lines indicate smoothed values fit by convolving the data 
with a Blackman window kernel filter. B) Numerical simulations under the non-dimensional model to fit X 1 278b (blue) and S288c (red). 
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Table 3 Parameter values used to replicate cAMP 
dynamics for strains £ 1 278b and S288c (both wild types) 
compared to parameter values used to fit the wild type 
data reported in Ma et al. [35] 



Parameter 


Interpretation 


Ma 


El 278b 


S288c 


M 


Rate at which Pde activity 
goes to steady state 


.01 


.085 


.19 


D 


Decay rate of cAMP due to 
active Pde1 


1 


0.247 


1.07 


r 


Pde1 affinity fore AMP 


33.6 


0.32 


0.15 


T 


Time scale 


.0377 


.17 


.35 



parameters for the pdelA and pde2A mutants we relaxed 
Condition (b), meaning that we no longer assumed a 
uniformly slower rate of activation/inactivation for Pde2. 
Model parameter estimates for the £ 1278b PDE mutants 
are given in Table 4. Comparing the parameters of the 
pdelA mutants in W303 (Table 1) and £ 1278b (Table 4), 
we see that the non-dimensional parameter M, represent- 
ing the rate at which phosphodiesterase activity reaches 
its steady state, increases relative to the wild-type in the 
£ 1278b background, reflecting the lack of a transient 
cAMP response in this strain. The difference between the 
backgrounds in terms of steady state cAMP levels follow- 
ing stimulus is primarily reflected in the relative change in 
the parameter D. 

Long-term dynamics 

If we consider our model over time scales longer than 
five minutes, our model the cAMP-PKA pathway predicts 



that cAMP levels may either go to a unique steady state 
or experience sustained oscillations. Garmendia-Torres et 
al. [42] predicted sustained oscillations in the concentra- 
tion of cAMP at intermediate stress levels. Since stress 
affects our model through parameter A (the activation rate 
of Ras-GTP), we choose A as our bifurcation parameter. 
At intermediate stress levels, for example A = .014, the 
long-term concentration of cAMP experiences sustained 
oscillations, as seen by the red curve in Figure 8. When 
stress levels are low, for example A = 1.4, the long-term 
concentration of cAMP remains at a high steady-state 
value, as seen by the blue line in Figure 8. When stress 
levels are high, for example A = .005, the long-term con- 
centration of cAMP remains at a low steady-state value, as 
seen by the black line in Figure 8. 

An intuitive biological explanation of this behavior 
is that at intermediate stress levels Ras (r) alternates 
between active (GTP-bound, r ss ^ 1) and inactive (GDP- 
bound, r ss ^ 0) states. This oscillation in Ras states causes 
oscillations in the concentration of cAMP. When stress 
levels are low, for example in the glucose fed state, Ras is 
almost completely in the Ras-GTP state, as in Ma et al. [35] 
Cases 1-4, forcing the long-term concentration of cAMP 
to steady state. Similarly, when stress levels are high, Ras 
is primarily in the Ras -GDP state, as in Ma et al. [35] Case 
5, forcing the long-term concentration of cAMP to a low 
steady state. 

Analysis of sustained oscillations 

Analysis of our model suggests that sustained oscillations 
are the result of PKA feedback through the Ras pathway 
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Figure 7 The signaling behavior of PDE mutants in the £ 1 278b strain background. A) Experimental time series of cAMP concentration 
following glucose stimulus for £ 1278b wild-type (blue), pdelA (red) and pde2A (green) strains, x's indicate measured cAMP values, the solid lines 
indicate smoothed values, as in Figure 7. B) Numerical simulations under the non-dimensional model to fit to the observed data. 
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Table 4 Parameter values used to replicate cAMP 
dynamics for the phosphodiesterase mutants, pdel A and 
pde2A in the £ 1278b background 



Parameter 




pc/el A 


pdel A 


M 


0.085 


0.262 


0.121 


D 


0.247 


0.0112 


0.107 


r 


0.32 


0.04 


0.061 



For all three cases, the time scale parameter T = 0.1 7. 



(variables r and z). In contrast, negative feedback through 
the PDEs (variable p) stabilizes the system. This analysis 
was done by examining the stability of the system by eval- 
uating the Jacobian matrix at the equilibrium. We find that 
our system loses stability through a Hopf bifurcation; thus, 
indicating that the system experiences sustained oscil- 
lations (for details see Additional file 1: Supplementary 
Model). This instability, leading to sustained oscillations, 
is a result of Ras varying between its active (GTP) and 
inactive (GDP) states. Activation/inactivation of Ras is 
controlled by parameters A, B> and Ti in our model. We 
find that as the value of Ti increases our system is less 
likely to experience sustained oscillations. Oscillations are 
also more likely as the ratio ^ increases. In regards to 
the period and amplitude of oscillations, we observe that 
as both forms of PKA feedback (through Ras and PDE) 
slow down, that is as the values of parameters N and M 
decrease, the length of the period of oscillation increases. 
Finally, we conclude that as the parameters that enhance 
the decay of cAMP are decreased, that is as parameters Do 



and D are decreased and parameter T is increased, oscil- 
lations are more likely to have longer periods and larger 
amplitudes. 

Discussion 

The model presented here is, to the best of our knowledge, 
the first analytical model of the yeast cAMP-PKA path- 
way that is capable of recapitulating both the short-term 
and long-term dynamics of cAMP signaling in wilde-type 
cells. The model we propose is also capable of replicat- 
ing cAMP dynamics that have been observed single and 
double phosphodiesterase mutants. In addition to cap- 
turing behaviors previously described our model makes 
new predictions about cAMP oscillations. Our work also 
highlights variation in cAMP signaling between genetic 
backgrounds. 

cAMP oscillations 

Analysis of a non-dimensional system derived from our 
model suggests that a substantial fraction of the model 
parameter space for the wild-type signaling network 
should result in decaying cAMP oscillations following 
glucose stimulus. While long time-scale oscillations have 
previously been predicted [39,42], short-term oscilla- 
tions, on the scale of minutes, have not been noted in 
the yeast literature. To test our predictions we experi- 
mentally measured the cAMP response in two different 
diploid strains of S. cerevisiae and observed oscillations 
as predicted by the model. While ours is the first study 
to document short time-scale cAMP oscillations in yeast, 
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Figure 8 Long-term dynamics for the concentration of cAMP. The blue line shows the cAMP concentration at high nutrient levels and/or low 

stress levels (parameter A = 1 .4). The black line shows the concentration of cAMP at low nutrient levels and/or high stress levels (parameter 

A = .005). The red curve shows the concentration of cAMP at intermediate nutrient levels and/or intermediate stress levels (parameter^ = .014). 
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such oscillations have been noted in other systems includ- 
ing mammals, amphibians, and slime molds [55-57]. For 
example, Dyachock et al [55] reported short time-scale 
glucose induced oscillations in mouse pancreatic /?-cells. 
They showed that these oscillations correlated with pul- 
satile insulin release. Similar oscillations are observed 
in response to stimulation of such cells by the peptide 
hormone GLP-1 [58]. The slime mold Dictyostelium also 
exhibits short time-scale cAMP oscillations [57]. Similar 
to the model we present here, Maeda et al. [57] suggested 
that positive feedback by PKA on phosphodiesterase 
activity and PKA negative feedback on adenylate cyclase 
activity are critical features of the molecular circuit that 
generates these oscillations in Dictyostelium, Whether 
cAMP oscillations in these diverse eukaryotic groups are 
due to conserved interactions in the cAMP pathway or 
whether such dynamics represent a convergent feature of 
cAMP signaling is a key question for future comparative 
studies. 

Our model is also able to generate longer-time scale 
cAMP oscillations. Such oscillatory dynamics should be 
reflected in PKA activity. Recent studies [59,60] have 
demonstrated that PKA activation plays an important role 
in regulating oscillatory patterns of nuclear localization 
of stress-responsive transcription factors such as Msn2. 
Different dynamical patterns of Msn2 regulation in turn 
can lead to qualitatively different expression outputs of 
downstream targets [59]. 

Comparison to other models 

There have been a number of recent studies that have 
proposed mathematical models for the cAMP-PKA net- 
work in yeast, including Gonze et al. [39], Cazzaniga et 
al. [40], and Williamson et al. [41]. A limitation shared by 
all three of these models is that they can not recapitulate 
the dynamical behaviors of the double phosphodiesterase 
mutant (pdelApde2A). This can be seen by setting the 
variables representing the phosphodiesterases to zero in 
any of these models (i.e. representing the double knock- 
out); the result is a rapid accumulation of cAMP. This 
contrary to the observed data reported by Ma et al. [35]. 
Our model uses competitive inhibition between Pdel, 
Pde2 and the GTPase activating proteins Iral and/or Ira2 
to explain the double PDE knockout case; thus, competi- 
tive inhibition between Pdel, Pde2, and Ira 1/2 is what sets 
our model apart from these previous models. 

Gonze et al. [39] propose a stochastic model, that is 
a direct extension of the Garmendia-Torres et al. model 
[42], to explore the effects of the cAMP-PKA path- 
way on the transcription factor Msn2. This model only 
focuses on the long term dynamics and can not be used 
to understand any of the PDE mutant cases proposed by 
Ma et al. Cazzaniga and colleagues [40] employed 
a stochastic modeling approach to study the yeast 



cAMP-PKA network. Though they do explore the effects 
of different PDE activities, their model does not explicitly 
account for the mutant cases described by Ma et al. and 
again they can not replicate the double phosphodiesterase 
knockout. Also they do not consider the long term 
dynamics proposed by Garmendia-Torres et al. [42]. The 
deterministic model presented by Williamson et al. [41] 
is the closest in approach to ours. Their model is able to 
recapitulate the relevant dynamics of the Pdel and Pde2 
single mutants but does not account for the behavior of 
the PDE double mutant. As we mentioned above this is 
significant since the double knockout case is the most 
surprising case and the hardest case to explain. They also 
only examined the short term dynamics of cAMP, where 
as our model accounts for both the long and short term 
dynamics of cAMP. 

Conclusions 

The mathematical model we present here replicates 
the dynamical behavior of wild-type and mutant cells 
and leads us to novel predictions of oscillatory cAMP 
dynamics in S. cerevisiae. Experiments motivated by this 
modeling effort confirmed these predictions and also 
highlighted variation in oscillatory phenomena among 
yeast isolates. The mathematical model, in turn, is capa- 
ble of reproducing the different signaling phenotypes we 
observed and allows us to predict which parts of the 
network are most likely to contribute to signaling vari- 
ation. Because it has both predictive and explanatory 
power this model can serve as a foundation for future 
mathematical and experimental explorations of cAMP 
signaling. A future challenge for both experimental and 
modeling studies of cAMP signaling will be to deter- 
mine whether the grossly similar oscillatory dynamics 
observed across eukaryotes is due to a common mech- 
anistic basis or whether this represents an example of 
convergent evolution. 

Methods 

Strains and cAMP assay 

Diploid strains of the laboratory backgrounds £ 1278b 
(G85 [wild- type', ura3<50, his350::hisG], G30025 [pdel A', 
pdel::bcKanMX4, ura3A0, his3A0::hisG], G30025 [pde2A\ 
pde2::bcKanMX4, ura3A0, his3A0::hisG], gifts of T. 
Galitski [61]) and S288c (BY4743, gift of D. Lew) were 
used to study cAMP dynamics. Cells were prepared 
following a protocol adapted from Paiardi et al. [17]. 
Approximately 2 x 10 9 cells/ml were grown in rich media 
(1% yeast extract, 2% peptone and 2% glucose) at 30°C 
and transferred to SC media with 0.1% glucose and 3% 
gycerol for overnight incubation at 30°C. Cells were col- 
lected, washed and incubated in 25 mM MES buffer, pH 
6 (Boston Bio Products, Worcester, MA) for 30 minutes 
at 30°C before glucose induction. Glucose was added to 
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2 ml of cell culture(approx. 2 xlO 8 cells/ml) to a final 
concentration of 100 mM and 25 mM. After glucose 
induction cells were collected over a 12-minute time 
course in 15 second increments for the first minute and 
at 30 second increments for the rest of the time course. 
5 x 10 6 cells were collected at each time point, fixed in 300 
/xl n-butanol-saturated 1M formic acid (Jiang etal, 1998) 
and frozen at -80°C. Cells were lysed using four freeze- 
thaw cycles and freeze-dried in a Speed Vac concentrator 
(Savant Instruments, Farmingdale, NY). Cell extracts 
were used to determine cAMP concentration using a 
cAMP Biotrak Enzyme immunoassay kit (Amersham, 
GE Healthcare) and a SpectraMax spectrophotometer 
(Molecular Devices, Sunnyvale, CA). 

Model fitting 

For the purposes of model fitting, where available we used 
estimates of key biochemical parameters as reported in 
[42]. In cases where such parameters were not avaiable, 
we used a method of least squares, implemented in the 
subroutine lmdif.f in the Fortran MINPACK library [62]. 
This subroutine minimizes the sum of the squares of m 
nonlinear functions in n variables by a modification of 
the Levenberg-Marquardt algorithm. We wrote a C com- 
puter program that calculates the appropriate functions 
that are used by this algorithm for the model described 
above. Additional details of the model fitting proce- 
dure are provided in the Additional file 1: Supplemental 
materials. 
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